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The Woods-Saxon basis has been suggested to replace the widely used harmonic oscillator basis for 
solving the relativistic mean field (RMF) theory in order to generalize it to study exotic nuclei. As 
examples, relativistic Hartree theory is solved for spherical nuclei in a Woods-Saxon basis obtained by 
solving either the Schrödinger equation or the Dirac equation (labelled as SRHSWS and SRHDWS, 
respectively and SRHWS for both). In SRHDWS, the negative levels in the Dirac Sea must be 
properly included. The basis in SRHDWS could be smaller than that in SRHSWS which will 
simplify the deformed problem. The results from SRHWS are compared in detail with those from 
solving the spherical relativistic Hartree theory in the harmonic oscillator basis (SRHHO) and those 
in the coordinate space (SRHR). All of these approaches give identical nuclear properties such as 
total binding energies and root mean square radii for stable nuclei. For exotic nuclei, e.g., ’Ca, 
SRHWS satisfactorily reproduces the neutron density distribution from SRHR, while SRHHO fails. 
It is shown that the Woods-Saxon basis can be extended to more complicated situations for exotic 
nuclei where both deformation and pairing have to be taken into account. 
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7 I. INTRODUCTION such as the Skyrme Hartree-Fock (SHF), Hartree-Fock- 
Bogoliubov (HFB), relativistic Hartree (RH) and rela- 


tivistic Hartree-Bogoliubov (RHB) theories. In these ap- 


j i proaches, solution of corresponding equations is trans- 
l vealed by the exceptional role of the nuclear magic num- formed to a matrix diagnolization problem which can be 
œ bers provides the foundation of the nuclear shell model 


— easily dealt with. 
€ and various mean field approaches [1, 2, 3]. This aver- 


age field is believed to be approximated most closely by 
a Woods-Saxon (WS) potential [4] either from analyzing 
the radial dependence of the nuclear force or by deriving 
it from a microscopic two-body force. 


The existence of an average field in atomic nuclei re- 


However, due to the incorrect asymptotic property of 
the HO wave functions, the expansion in the localized 
HO basis are not appropriate for the description of drip 
line nuclei [6, 7, 8] which display many interesting fea- 
tures because of the extremely weakly bound property, 


Since the eigenfunctions for the WS potential can not 
be given analytically, as good approximations for stable 
nuclei, one often adopts the harmonic oscillator (HO) 
potential or the square well, in particular the former, 
in shell model calculations for both spherical [1] and de- 
formed nuclei [5]. The HO eigenfunctions also often serve 
as a complete basis in solving equations in both non- 
relativistic and relativistic mean field approximations, 
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e.g., the coupling between bound states and the contin- 
uum due to the pairing correlation, large spacial density 
distributions, possible modifications of shell structure, et 
al. One must improve the asymptotic behavior of HO 
wave functions, e.g., by performing a local scaling trans- 
formation [9]. However, one does not know the scaling 
parameter beforehand thus the predictive power in this 
method is lost. 


A proper representation to solve the HFB or RHB 
equations for drip line nuclei is the coordinate space [7, 
10, 11, 12] where wave functions are approximated on a 
spatial lattice and the continuum is discretized by suit- 
ably large box boundary conditions. The HFB method 
solved in r space can take fully into account all the mean- 
field effects of the coupling to the continuum [6, 7, 10, 13]. 


Nevertheless for deformed nuclei, working in r space be- 
comes much more difficult and numerically very sophis- 
ticated [8]. Particularly, it become very time consuming 
when the pairing correlation is included. Therefore much 
effort is made towards a more efficient solution of HFB or 
RHB equations, e.g., using natural orbitals [14] or work- 
ing on basis-spline Galerkin lattices [15], etc.. 

A reconciler between the HO basis and the r space 
may be the WS basis because (i) the WS potential repre- 
sents the nuclear average field more suitably than the HO 
potential and (ii) in principal there is no localization re- 
strictions in the WS potential. Although analytical wave 
functions can not be given for the WS potential, one may 
easily find numerical solutions for a spherical WS poten- 
tial in the r space by virtue of various effective methods 
of solving ordinary differential equations [16]. One can 
still use a large box boundary condition to discretize the 
continuum. These WS wave functions can thus be used 
as a complete basis for spherical or deformed systems and 
one finally comes back to the familiar matrix diagnoliza- 
tion problem. 

In the present work we restrict the application of this 
method to nuclei with spherical symmetry which largely 
facilitates the discussion of basic principles and allows 
presenting illustrations for the radial dependence of all 
relevant physical quantities like density distributions. We 
combine this approach with the relativistic Hartree the- 
ory [17] which provides a framework for describing the 
nuclear many body problem as a relativistic system of 
baryons and mesons and, together with its extensions 
with deformation and/or pairing included, have been suc- 
cessfully applied in calculations of nuclear matter and 
properties of finite nuclei throughout the periodic ta- 
ble [18, 19]. 

The paper is organized as follows. In Sec. II, we give 
a brief reminder of the formalism of relativistic Hartree 
theory. The numerical details of solving it in the WS 
basis are given in Sec. III. In Sec. IV we present our 
results and compare them with those obtained in the HO 
basis and in the r space. We also discuss the contribution 
from negative levels in the Dirac sea in the same section. 
Finally, the work is summarized in Sec. V. 

Throughout the paper, the relativistic Hartree theo- 
ries solved in the r space, in the HO basis and in the 
WS basis are abbreviated as “SRHR”, “SRHHO” and 
“SRHWS” where the first “S” represents “spherical”. 
We use “SWS” and “DWS” to distinguish the WS ba- 
sis which is obtained from solving the Schrodinger equa- 
tion or the Dirac equation with initial WS potentials, 
respectively. Thus we have “SRHSWS” and “SRHDWS” 
theories. 


II. BASIC FORMALISM OF RELATIVISTIC 
HARTREE THEORY 


The starting point of the relativistic Hartree theory 
is a Lagrangian density where nucleons are described as 


Dirac spinors which interact via the exchange of several 
mesons (a, w, and p) and the photon [17, 18, 19], 


£ = Vili — M) hi + 50,00"0 — U(0) — goelio 
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with the summation convention used and the summation 
over i runs over all nucleons, £ = y £, = Y £“, M the 
nucleon mass, and Mo, go, Mw, Ju; Mp, Jp Masses and 
coupling constants of the respective mesons. The non- 
linear self-coupling for the scalar mesons is given by [20] 


1 
U(o) = zoo + r + Sot, (2) 
and field tensors for the vector mesons and the photon 


fields are defined as 
Quy = pwy — Irwy, 


Rw = Ou Pv = Ov Pu = Ip (Pu x pv); (3) 
Fav = OpAy — OA. 


The classical variation principle gives the equations of 
motion for the nucleons, mesons and the photon. As in 
many applications, we study the ground state properties 
of nuclei with time reversal symmetry, thus the nucleon 
spinors are the eigenvectors of the stationary Dirac equa- 
tion 


[a-p+V(r) + 6(M + S(r))] yir) = etir), (4 


and equations of motion for the mesons and the photon 
are 


(-A+0,U(c))o(r) = —gops(r), 
-A +m?) w? (r) = Ju P(r), (5) 
(-A +m?) P(r) = gpps(r), 
~AA"r) = epp(n). 


where w? and A? are time-like components of the vec- 
tor w and the photon fields and p° the 3-component of 
the time-like component of the iso-vector vector p me- 
son. Equations (4) and (5) are coupled by the vector and 
scalar potentials 
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2A (r), 


[vo = guwP(r) + gptsp%(r) +e (6) 


S(r) = gso(r), 


and various densities 


ps(r) = Di Bilr)bi(r), 

polr) = Dia ti lr)yilr), 

p3(r) = Zi Wr etlr), (7) 
pelr) = TA, dh (r) —Bui(r). 


For spherical nuclei, meson fields and densities depend 
only on the radial coordinate r. The spinor is character- 
ized by the angular momentum quantum numbers (l,j), 


m, the parity, the isospin t = +1/2 (“+” for neutrons 
and “—” for protons), and the radial quantum number 
a. The Dirac spinor has the form 
Et) v0, ¢) 1 
Panm(r,s,t) = FF(r) x Xta (t), J =lt5, 
a y}n(8, 4) 


(8) 
with Gk(r)/r and F¥(r)/r the radial wave functions for 
the upper and lower components and Ys (0, p) the spin 
spherical harmonics where « = (—1)/+'+1/2(j + 1/2) and 
i =1+(-1)9+'-1/2_ The value of « of the upper compo- 
nent is used to label a state both for normal levels in the 
Fermi sea and for negative ones in the Dirac sea. States 
with the same « form a “block”. The radial equation of 
the Dirac spinor, Eq. (4), is reduced as 


€aGG = Or r FE + (M + 8S(r)+V(r)) GS, 
ee = ee NAME E 
Or r 


(9) 
The meson field equations become simply radial Laplace 
equations of the form 


8? 
a 
Mg are the meson masses for ¢ = o,w, p and zero for the 
photon. The source terms are 


2am’) ofr) = 0) 


—~gops(r) — g20?° (r) — gao? (r), for a, 
Ju Pv (r), for W, 


907) = 9 gppslr), for p, OH) 
epe(r), for A, 
with 
Ane pdr) = Lile) - |F;(r)|?), 
4nr?py(r) = Zim (lG: (r)? + |F;(r)|?), (12) 
Anr?p3(r) = diet 2t;(|Gi(r)|? + |Fi(r)|?), 
4rr?pelr) = Yh- (3 — ti) (Gs)? + |B). 


The above coupled equations have been solved in the r 
space [21] and in the HO basis [22] with the no sea and the 
mean field approximation. Here we depict briefly the pro- 
cedure of solving these coupled equations. With a set of 
estimated meson and photon fields, the scalar and vector 
potentials are calculated and the radial Dirac equation is 
solved. Thus obtained nucleon wave functions are used 
to calculate the source term of each radial Laplace equa- 
tion for mesons and the photon. New meson and photon 
fields are calculated from solving these Laplace equations. 
This procedure is iterated until a demanded accuracy 
is achieved. Laplace equations are usually solved using 


TABLE I: Dependence of the average single particle energy, 
root mean square (rms) radius, and (r*)'/* on the mesh size 
Ar for the SRHSWS theory. The meson and Coulomb fields 
are obtained from SRHR calculations with the parameter set 
NL3, Ar = 0.05 fm, Rmax = 30 fm for 16O and 35 fm for 
208Pb. In SRHSWS calculations, the parameter set NL3 is 
used. For '©O, Rmax = 4roA!/? = 12.8 fm, Ecut = 300 MeV. 
For 2°°Pb, Rmax = 3ro Al? = 22.6 fm, Ecut = 200 MeV. The 
first row gives results from solving the Dirac equation in the 
coordinate space. 


Ar —Esp/A (21/2 (rt) 1/4 


(r4) 1/4 


O Pb 
— | 23.0375 2.5945 2.8900] 23.3348 5.6315 5.9883 
0.05} 23.0375 2.5945 2.8900] 23.3348 5.6315 5.9883 
0.10} 23.0376 2.5945 2.8899] 23.3348 5.6315 5.9883 
0.20} 23.0385 2.5941 2.8896] 23.3347 5.6315 5.9883 
0.30} 23.0430 2.5930 2.8889] 23.3343 5.6314 5.9883 
0.40} 23.0420 2.5923 2.8885] 23.3334 5.6315 5.9885 
0.50} 22.9887 2.5949 2.8936] 23.3284 5.6319 5.9890 


the Green’s function method [21, 22] though in Ref. [22] 
Laplace equations for mesons are solved in the HO ba- 
sis. SRHR, SRHHO and SRHWS differ from each other 
mainly in how to solve the Dirac equation. In the follow- 
ing the numerical solution of the Dirac equation in the 
WS basis will be presented. 


III. SOLVING THE DIRAC EQUATION IN A 
WOODS-SAXON BASIS AND NUMERICAL 
DETAILS 


A. Woods-Saxon basis from solving a Schrodinger 
equation (the SWS basis) 


For the Schrodinger equation with a spherical Woods- 
Saxon potential 


Vo 


Vws(r) = 1 + e(r—Ro)/a0’ 
OO, r > Rmax; 


r < Rmax, (13) 


where Rmax is introduced for practical reasons to define 
the box boundary, the eigenfunction can be written as 
nim (Tr) = Rni(r)Yim, (0, $). Its radial Schrödinger equa- 
tion is derived as 


-ar (SRF H) Hio] rae) 
= EmRnlr). (14) 


Equation (14) is solved on a discretized radial mesh 
with a mesh size Ar. Rmax (Ar) should be chosen larger 
(smaller) enough to make sure that the final results do 
not depend on it. The radial wave functions thus ob- 
tained form a complete basis, 


{Rni(r);n =0,1,--- ;1=0,1,--- ,n}, (15) 


TABLE II: Dependence of the average single particle energy, 
rms radius, and (43/4 on the difference An = ftmax — Nmax 
for the SRHSWS theory. The meson and Coulomb fields 
are obtained from SRHR calculations with the parameter set 
NL3, Ar = 0.05 fm, Rmax = 30 fm for 16O and 35 fm for 
208Pb. In SRHSWS calculations, the parameter set NL3 is 
used. For 160) Rast = Arp A13 = 12.8 fm. For 208Pp, 


Rmax = 3ro A13 = 22.6 fm. 
An —Esp/A (r?) 1/2 (r4) 1/4 _ sp/A (r?) 1/2 (r4) 1/4 
Pb: Ecut = 100 MeV 


O: Ecut = 100 MeV 


1 | 23.0382 2.5947 2.8920] 23.3344 5.6315 5.9884 
3 | 23.0326 2.5949 2.8913] 23.3341 5.6315 5.9884 
5 | 23.0298 2.5952 2.8915] 23.3340 5.6315 5.9884 
7 | 23.0290 2.5953 2.8915] 23.3340 5.6315 5.9884 
9 | 23.0289 2.5953 2.8915] 23.3340 5.6315 5.9884 

O: Ecut = 300 MeV Pb: Ecut = 200 MeV 
1 | 23.0375 2.5945 2.8999] 23.3348 5.6315 5.9883 
3 | 23.0375 2.5945 2.8999] 23.3348 5.6315 5.9883 
5 | 23.0375 2.5945 2.8999] 23.3347 5.6315 5.9883 
7 | 23.0375 2.5945 2.8999] 23.3347 5.6315 5.9883 
9 | 23.0375 2.5945 2.8999] 23.3347 5.6315 5.9883 


in terms of which the radial part of the upper and the 
lower components of the Dirac spinor in Eq. (9) are ex- 
panded respectively as 


G(r) =—i 5 Jant Rni(r), 
n=0 


Nmax 


FEl) =i X foarRgi(r). 


n=0 


(16) 


The radial Dirac equation, Eq. (9), is transformed into 
the WS basis as 


Amn Bra Jan no Jan 
( Crn Dinn ) ( Jañ ) = €a ce ) ; (17) 


where the matrix elements are calculated as follows 


(18) 
In practical calculations, an energy cutoff Eeut (relative 
to the nucleon mass M) is used to determine the cutoff 
of the radial quantum number nmax for each block. In 
the expansion of the corresponding lower component, we 
take fmax = Nmax + An with An > 1 in order to avoid 
spurious states [22]. 
The following Woods-Saxon parameters have been 
used according to Ref. [24] 


Vo = (—51 + 33(N — Z)/A) MeV, 
Ro = 1.27A'/3 fm, ag = 0.67 fm, 


(19) 
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FIG. 1: Deviations of the average single particle energy 


E*?/A (the upper panel), rms radius (r)!/? (the middle 
panel), and (r4)!/4 (the lower panel) of 16O from the standard 
results versus the cutoff energy Ecut with different box size 
Rmax for the SRHSWS theory. The meson and Coulomb fields 
are obtained from SRHR calculations with the parameter set 
NL3, Ar = 0.05 fm, Rmax = 30 fm. In SRHSWS calculations, 
the parameter set NL3 is used. For Ecut = 100, 200, 300, 400 
MeV, the approximate maximum principal quantum number 
in each basis, Nmax = 2Mmax + l, is given on the top of the 
plot. 


where ‘+’ is for the neutron and ‘—’ for the proton. As ex- 
pected, the dependence of final results on the initial WS 
potential is almost negligible. For example, a variation 
of Vo by 50% gives differences in total binding energies 
by less than 0.1% and charge radii by less than 0.5% for 
160, 48Ca and 2°8Pb. Such situation is also checked to be 
true for the other two parameters in the WS potential, 
Ro and ao. 


B. Woods-Saxon basis from solving a Dirac 
equation (the DWS basis) 


The radial Dirac equation, Eq. (9), may be solved in 
the r space [21] with Woods-Saxon-like potentials for 
Vo(r)+So(r) [23] within a spherical box of the size Rmax, 
together with the spherical spinor which gives a complete 
WS basis 


{ Ennis T (r, 5, t)] ; ee = 0} ’ (20) 
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FIG. 2: Deviations of the average single particle energy 


E‘?/A (the upper panel), rms radius (r)'/? (the middle 
panel), and (r*)'/4 (the lower panel) of 7°°Pb from the stan- 
dard results versus the cutoff energy Ecut with different box 
size Rmax for the SRHSWS theory. The meson and Coulomb 
fields are obtained from SRHR calculations with the parame- 
ter set NL3, Ar = 0.05 fm, Rmax = 35 fm. In SRHSWS calcu- 
lations, the parameter set NL3 is used. For Ecut = 100, 200, 
300, 400 MeV, the approximate maximum principal quantum 
number in each basis, Nmax = 2Nmax + l, is given on the top 
of the plot. 


with n = 0,1,---,® = +1,+2,---, and M = —jxy- ++ Ju. 

O mlr, 5,t) takes the form of Eg: (8). We note that 
states both in the Fermi sea and in the Dirac sea should 
be included in the basis for the completeness. The nu- 
cleon wave function, Eq. (8), can be expanded in terms 


of this set of basis as 


Nmax 


= Cog elt, S st); (21) 


Warm T, s,t) 


where Nmax = nax nmax +1 and the seen is over 
normal levels in the Fermi sea for 0 < n < nif, and over 
negative levels in the Dirac sea for niay +1 <n < Nmax: 
The negative states is obtained with the same method 
as the positive ones [21]. In this WS basis, the Dirac 
equation, Eq. (4), turns out to be 


Nmax 


f 
Camé?, + > Can H. 


n=0 


> Nmax, (22) 
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FIG. 3: Deviations of the average single particle energy 


E®/A (the upper panel), rms radius (r)!/? (the middle 
panel), and (r*)!/4 (the lower panel) of 16O from the stan- 
dard results versus the cutoff energy Ecut with different box 
size Rmax for the SRHDWS theory. The meson and Coulomb 
fields are obtained from SRHR calculations with the parame- 
ter set NL3, Ar = 0.05 fm, Rmax = 30 fm. In SRHDWS cal- 
culations, the parameter set NL3 is used. For Ecut = 100, 200, 
300, 400 MeV, the approximate maximum principal quantum 
number in each basis, Nmax = 2Nmax + L, is given on the top 
of the plot. 


with 
Hin = (Ymr) | [AV (r) + BAS(r)] [YO (r 
= | ion dr@®, (r) [AV (r) + AS(r)] 7 (r) 
0 


Rass 
+f dr F® (r) [AV (r) — AS(r)] F? (r) (23) 
0 


where AV (r) = V(r)—Vo(r) and AS(r) = S(r)— So(r). 
The angular, spin, and isospin quantum numbers are 
omitted for brevity. 

It should be mentioned that Eq. (9) can be solved di- 
rectly in the r space with the same method of generating 
the DWS basis. It is our aim to test the validity of an 
efficient solution not only for the spherical RH model but 
also for its extension to include the deformation and/or 
the pairing correlation. In fact, if only the SRH theory 
is concerned, this procedure is just a replacement of the 
direct solution in the r space by a diagnolization of a 
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FIG. 4: Deviations of the average single particle energy 


E‘?/A (the upper panel), rms radius (r)'/? (the middle 
panel), and (r*)'/4 (the lower panel) of 7°°Pb from the stan- 
dard results versus the cutoff energy Ecut with different box 
size Rmax for the SRHDWS theory. The meson and Coulomb 
fields are obtained from SRHR calculations with the parame- 
ter set NL3, Ar = 0.05 fm, Rmax = 35 fm. In SRHDWS cal- 
culations, the parameter set NL3 is used. For Ecut = 100, 200, 
300, 400 MeV, the approximate maximum principal quantum 
number in each basis, Nmax = 2Nmax + l, is given on the top 
of the plot. 


matrix with some complication introduced by the fact 
that contributions from states in the Dirac sea must be 
included. 

An energy cutoff Ecut (relative to the nucleon mass 
M) and the cutoff of radial quantum numbers nj, are 
applied to normal levels alternatively according to prac- 
tical convenience. For the initial Woods-Saxon potentials 


Vo(r) + So(r), we follow Ref. [23]. 


C. Comparison with the r-space method 


In order to check the validity of solving the Dirac equa- 
tion in the WS basis and to provide numerical experiences 
for future applications, we compare results of 160 and 
208Pb from solving the Dirac equation in the WS basis 
and those from solving the same equation in the coordi- 
nate space. The latter is the most accurate method of 
solving the Dirac equation for realistic nuclei up to now 


thus is used as a standard here. The scalar and vector 
potentials in the Dirac equation are provided by very ac- 
curate SRHR calculations with the parameter set NL3 
for the Lagrangian, the mesh size Ar = 0.05 fm, the box 
size Rmax = 30 fm for 160 and Rmax = 35 fm for SPb. 
Then with thus obtained S(r) and V(r), the Dirac equa- 
tion is solved in the coordinate space and in the WS basis 
with also the parameter set NL3. 

To compare results only from solving the Dirac equa- 
tion avoids errors from other numerical procedures, e.g., 
the error from the iteration and that from solving the 
Laplace equations. For the same reason, what we com- 
pare between these two methods is not the binding energy 
which contains the contribution from mesons but the av- 
erage single particle energy E*P/A. ESP = 5°, e; where e; 
is the single particle energy and the summation over all 
occupied states for both neutrons and protons. We also 
compare the rms radius (r?)!/? and (r*)1/4. The radius 
(r4)1/4 reflects the nucleon densities in the tail more than 
the rms radius. 

Table I presents the dependence of results of the Dirac 
equation in the SWS basis on the mesh size Ar. With 
Ar decreasing, results in the SWS basis approach the 
standard results, i.e., those in the r space. Ar = 0.1 fm 
gives results accurate enough. The dependence on the 
box size Rmax and on the basis size determined by Ecut 
are investigated and shown in Figs. 1 and 2 where the 
deviations of the average single particle energy E*?/A, 
the rms radius (r?)!/? and (r+)1/4 from the standards 
are plotted versus Ecut for different Rmax. If Rmax is 
not large enough, it is difficult to approach the standard 
results. For example, when Rmax = 3roA!/? = 9.4 fm 
and Eeut = 300 MeV for 160 (correspondingly, Nmax ~ 
20), the results seem converge, but the discrepancy of 
the average single particle energy from the standard one 
remains 0.1 keV. So one must use a large enough box with 
the size Rmax around 479 A'/ for light nuclei and 3r9 A!/? 
for heavy ones. It is interesting that the convergence of 
the results does not depend on Nmax but only on Ecut. 
For 160 (?°8Pb), the results converge to the standard 
ones at ~ 300 (400) MeV. From Figs. 1 and 2, we find 
the radius (r*)!/4 also converges very well which implies 
that nucleon densities could be calculated accurately even 
for large r. 

We have made similar investigations for results in the 
DWS basis and similar conclusions are made. For in- 
stances, the deviations of the average single particle en- 
ergy ESP /A, the rms radius (r?)!/? and (r+)1/4 from the 
standards are plotted versus Eecut for different Rmax in 
Figs. 3 and 4. 

In the expansion of the nucleon wave function, 
Eq. (21), one has to take into account not only the lev- 
els in the Fermi sea but also those in the Dirac sea be- 
cause they form a complete basis together. Now the ques- 
tion arises how many levels in the Dirac sea one has to 
take into account. In the calculations in Figs. 3 and 4, 
we have used nias = nhay With ntas determined by 
FEeut- In Table III, the dependence of the average sin- 


TABLE III: Dependence of the single particle energy, rms ra- 
dius, and (r41 t on the maximum principal quantum number, 
Nmnax = 2Mnax + l, for the SRHDWS theory. The meson and 
Coulomb fields are obtained from SRHR calculations with the 
parameter set NL3, Ar = 0.05 fm, Rmax = 30 fm for 1°O and 
35 fm for 2°°Pb. In SRHDWS calculations, the parameter set 
NL3 is used. For 1®O, Rmax = 4r0A'/? = 12.5 fm, Ecut = 300 
MeV for positive states. For 208 Pp, Rae = 3ro A13 = 22.1 
fm, Ecut = 200 MeV for positive states. 


Nmax| -Bool A (YP (A ] -EA ye ey” 
no | 23.1129 2.5912 2.8859] 23.3331 5.6314 5.9889 
0 23.1077 2.5916 2.8861] 23.3329 5.6314 5.9889 
2 23.0762 2.5939 2.8889] 23.3316 5.6315 5.9890 
4 23.0617 2.5942 2.8893] 23.3304 5.6317 5.9892 
6 23.0439 2.5946 2.8898 | 23.3299 5.6318 5.9893 
8 23.0385 2.5946 2.8899] 23.3294 5.6319 5.9893 
10 | 23.0376 2.5946 2.8899] 23.3292 5.6319 5.9894 
12 | 23.0375 2.5946 2.8899] 23.3291 5.6319 5.9899 
14 | 23.0375 2.5946 2.8899] 23.3290 5.6319 5.9894 
16 | 23.0375 2.5946 2.8899] 23.3290 5.6319 5.9894 
18 | 23.0375 2.5946 2.8899] 23.3289 5.6319 5.9893 
20 | 23.0375 2.5946 2.8899} 23.3288 5.6319 5.9893 
22 | 23.0375 2.5946 2.8899} 23.3287 5.6319 5.9893 
30 | 23.0375 2.5946 2.8899} 23.3287 5.6319 5.9893 


gle particle energy, the rms radius (r?)!/? and (r*)!/4 on 
N nax = 2%max +! — a cutoff on the principal quantum 
number of levels in the Dirac sea — are given for 160 
and 2°8Pb, From Table III, we find a merit of solving the 
Dirac equation in the DWS basis: the number of nega- 
tive states included in the basis could be much smaller 
that that of the positive states. Let take 16O as an ex- 
ample, Rmax = 4ro A13 and Ecut = 300 MeV for positive 
states correspond to Nt, ~ 28. But for negative states, 
Nwax = 10 gives very accurate results, e.g., the discrep- 
ancy of ESP /A from the standard is smaller than 0.1 keV. 
This will significantly simplify the deformed problem by 
decreasing the matrix dimension compared to solve the 
Dirac equation in the SWS basis. 

The above investigations are somehow academic. In 
practical applications, it is not necessary to go to the 
accuracy around keV in the single particle energy or 1074 
fm in the radius. So in the following calculations, we 
will use Rmax = 20 fm, Ar = 0.1 fm and Egy; = 60~80 
MeV for heavy and light nuclei which give reasonable 
accuracies both for the binding energy and the radius. 
This set of cutoff’s corresponds approximately Nmax = 
2Nmax +l ~ 25 where | is the orbital angular momentum 
of relevant state. 


IV. RESULTS AND DISCUSSIONS 


In this section we present results of SRHWS. Since 
our main aim is to show the virtues of SRHWS com- 
pared to SRHHO and SRHR, we do not include pairing 


correlations and restrict our study to doubly magic or 
magic nuclei only. If not specified, the parameter set 
NLSH is used for the Lagrangian, Rmax = 20 fm and 
Ar = 0.1 fm throughout this section. Other parame- 
ter sets for the Lagrangian do not change the conclu- 
sion here. In SRHDWS, the number of normal levels in 
the Fermi sea and that of negative ones in the Dirac sea 
are the same for convenience, i.e., Nia, = Nmax For 
SRHHO, hwy = 41A~!/3 has been used and cutoff’s of 
the expansion for fermions and bosons are the same, i.e., 
Ne = Ng = Nane 


A. Bulk properties of stable nuclei from different 
SRH theories 


In Table IV, the binding energy per nucleon (E/A) and 
neutron, proton and charge radii (rn, rp and re) of some 
typical spherical nuclei are presented which are calcu- 
lated from the present available codes, including SRHR, 
SRHSWS, SRHDWS, SRHHO. Available data [25, 26] 
are also included for comparison. We use approximately 
the same Nmax in the SRHHO as that in the SRHWS 
which is determined by Eecut- 


Generally speaking, for each studied nucleus, the four 
approaches give almost the same results with an accuracy 
within 0.1% with few exceptions where the differences are 
still less than 0.3%. They are in excellent agreement with 
available data. 


With the same parameters of spatial lattice Rmax and 
Ar, SRHWS should reproduce results of SRHR when 
Ecut (or Nmax) is large enough. This is justified in 
Table IV. One find exactly coincident results between 
SRHSWS and SRHR for most of the studied nuclei. The 
remaining differences and those between SRHDWS and 
SRHR could be diminished by increasing Eeut. 


In Figures 5, 6 and 7, the neutron density distributions 
are compared between SRHR, SRHSWS, SRHDWS, and 
SRHHO, in which 160, #8Ca and ?°°Pb are chosen as 
examples. The calculation details are the same as Ta- 
ble IV. For these stable nuclei, all these SRH methods 
are valid and all calculations are in excellent agreement 
with each other from the central to outer region of each 
nucleus. Tiny differences in the central region do not af- 
fect the physical observables like the binding energy or 
nuclear radius as is seen in Table IV. Furthermore, these 
differences could also be decreased by increasing Ecut or 
Nmax 


From the above discussions, it is clear that SRHWS 
is equivalent to SRHR and SRHHO for stable nuclei. 
Thus we conclude that Woods-Saxon basis provide an- 
other possibility to solve (non-)relativistic mean field the- 
ory. 


TABLE IV: The binding energy per nucleon and neutron, 
proton and charge radii of some typical spherical nuclei. The 
parameter set NLSH is used for the Lagrangian. Rmax = 20 
fm and Ar = 0.1 fm for SRHR and SRHWS. hwo = 4147"? 
for SRHHO. Numbers in brackets in the second column give 
Ecut for SRHWS and Nmax for SRHHO. Data for E/A and 
re are taken from Ref. [25] and Ref. [26], respectively. Energy 
is in MeV and radius in fm. 


Nucleus E/A Tn Tp Te 

TO SRHR —8.022 2.551 2.578 2.699 
SRHSWS (80) -—8.022 2.554 2.581 2.702 
SRHDWS (80) -—8.014 2.553 2.580 2.701 
SRHHO (25) —8.034 2.551 2.577 2.699 
Experiment —7.976 2.693 

Ca SRHR —8.500 3.311 3.359 3.452 
SRHSWS (80) -—8.499 3.310 3.358 3.452 
SRHDWS (80) -—8.497 3.312 3.359 3.453 
SRHHO (25) —8.514 3.3810 3.358 3.452 
Experiment —8.551 3.478 

Ca SRHR —8.644 3.586 3.369 3.463 
SRHSWS (80) —8.646 3.583 3.368 3.461 
SRHDWS (80) -—8.639 3.586 3.371 3.464 
SRHHO (25) —8.659 3.584 3.368 3.462 
Experiment —8.666 3.479 

56Ni SRHR —8.634 3.582 3.630 3.717 
SRHSWS (80) -—8.640 3.580 3.628 3.715 
SRHDWS (80) -—8.625 3.585 2.633 3.720 
SRHHO (25) —8.650 3.581 3.629 3.716 
Experiment —8.345 

Zr SRHR —8.677 4.294 4.186 4.262 
SRHSWS (75) —8.677 4.295 4.187 4.263 
SRHDWS (75) —8.672 4.295 4.187 4.262 
SRHHO (25) —8.693 4.293 4.185 4.261 
Experiment —8.710 4.270 

TSSn SRHR —8.466 4.743 4.553 4.623 
SRHSWS (70) —8.466 4.743 4.554 4.624 
SRHDWS (70) —8.460 4.743 4.554 4.624 
SRHHO (25) —8.482 4.741 4.552 4.622 
Experiment —8.517 4.641 

T3ZSn SRHR —8.377 4.964 4.636 4.704 
SRHSWS (70) -—8.377 4.964 4.637 4.704 
SRHDWS (70) —8.370 4.964 4.637 4.706 
SRHHO (25) —8.393 4.963 4.635 4.703 
Experiment —8.355 

20Pp SRHR —7.885 5.713 5.447 5.505 
SRHSWS (60) —7.886 5.712 5.447 5.505 
SRHDWS (60) —7.874 5.712 5.448 5.506 
SRHHO (25) —7.900 5.711 5.445 5.504 
Experiment —7.868 5.504 


B. Neutron density distributions for “Ca in 
different SRH theories 


As we already discussed in the introduction, one of the 
merits of SRHR against SRHHO is its proper description 
of exotic nuclei. In this subsection, we will demonstrate 
the equivalence between SRHWS and SRHR when rea- 
sonably large Eecut is applied in SRHWS. 


In order to see the results for the unstable nuclei near 


0.10 


0.08 


— 0.06 
£ 
€ 
Æ 0.04 
— SRHR 
_ - SRHDWS 
0.02 ——- SRHSWS 
—-— SRHHO 
0.00 
0 3 6 
r [fm] 


FIG. 5: Neutron density distributions for 16O from different 
SRH approaches. The parameter set NLSH is used for the 
Lagrangian. Rmax = 20 fm and Ar = 0.1 fm for SRHR and 
SRHWS. Ecut = 80 MeV for SRHWS. Correspondingly, Nmax 
= 25 for SRHHO. In SRHDWS, the number of levels in the 
Dirac sea included in each block is the same as that of nor- 
mal levels which is determined by Ecut. The inset presents 
logarithmic densities. 
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FIG. 6: Neutron density distributions for “*Ca from differ- 
ent SRH approaches. The parameter set NLSH is used for 
the Lagrangian. Rmax = 20 fm and Ar = 0.1 fm for SRHR 
and SRHWS. Ecut = 80 MeV for SRHWS. Correspondingly, 
Nmax = 25 for SRHHO. In SRHDWS, the number of levels 
in the Dirac sea included in each block is the same as that of 
normal levels which is determined by Fcut. The inset presents 
logarithmic densities. 


the neutron drip line, the neutron density distribution for 
Ca is studied here. The nucleus ‘?Ca is predicted to 
be the last bound calcium isotope [27, 28, 29, 30]. Since 
it is not a doubly magic nucleus, there might be some 
uncertainty in present results due to the lack of inclusion 
of pairing correlations. However, as we stressed in the 
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FIG. 7: Neutron density distributions for ?°°Pb from differ- 
ent SRH approaches. The parameter set NLSH is used for 
the Lagrangian. Rmax = 20 fm and Ar = 0.1 fm for SRHR 
and SRHWS. Ecut = 60 MeV for SRHWS. Correspondingly, 
Nmax = 25 for SRHHO. In SRHDWS, the number of levels 
in the Dirac sea included in each block is the same as that of 
normal levels which is determined by FEcut. The inset presents 
logarithmic densities. 


TABLE V: Convergence study for Ca. The parameter set 
NLSH is used for the Lagrangian. Ar = 0.1 fm for SRHR and 
SRHWS. Ecut = 75 MeV for SRHWS. Energy is in MeV and 
radius in fm. 


E/A Tn Tp An 

Rmax SRHR 

20 6.482 4.656 3.639 —0.191 

25 6.483 4.723 3.639 —0.221 

30 6.484 4.773 3.639 —0.228 

35 6.484 4.807 3.639 —0.229 
Rmax SRHSWS 

20 6.481 4.663 3.639 —0.206 

25 6.482 4.726 3.639 —0.231 

30 6.483 4.774 3.639 —0.237 

35 6.483 4.803 3.639 —0.238 
Rmax SRHDWS 

20 6.474 4.662 3.641 —0.163 

25 6.475 4.733 3.641 —0.197 

30 6.475 4.789 3.640 —0.205 

35 6.475 4.828 3.640 —0.206 
Nmax SRHHO 

25 6.489 4.577 3.639 —0.054 

31 6.492 4.605 3.639 —0.128 

37 6.494 4.628 3.639 —0.166 

43 6.494 4.649 3.639 —0.189 


beginning of this section, the main aim here is to show 
the virtue of SRHWS compared to SRHHO, it is very 
unlikely that the pairing would change our conclusion 
qualitatively. 

For stable nuclei, it has been shown that Rmax ~ 20 
fm is large enough. For drip line nuclei, the dependence 
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FIG. 8: Comparison of density distributions for Ca from 
SRHR, SRHSWS and SRHHO. The parameter set NLSH 
is used for the Lagrangian. Ar = 0.1 fm for SRHR and 
SRHSWS. Ecut = 75 MeV and Rmax = 20 (thick dot-dashed 
curve), 25 (thick long-dashed curve), 30 (thick dashed curve) 
and 35 fm (thick dotted curve) for SRHSWS. These sets of 
cutoff’s correspond to cutoff’s in principal quantum number 
Nmax = 25 (thin dot-dashed curve), 31 (thin long-dashed 
curve), 37 (thin dashed curve) and 43 (thin dotted curve) 
which are used in SRHHO calculations. The density distribu- 
tion from SRHR are almost identical with that from SRHSWS 
with the same box size. For brevity, only pn(r) from SRHR 
with Rmax = 35 fm (thick solid line) is displayed here which 
covers the curve corresponding to pn(r) from SRHSWS with 
Rmax = 35 fm (thick dotted curve). 


of the results on Rmax for ‘*Ca is presented in Table V. 
For both SRHR and SRHWS, Ar = 0.1 fm and Rmax 
= 20, 25, 30 and 35 fm have been adopted respectively. 
The energy cutoff Eeut = 75 MeV is used to SRHWS. In 
SRHR and SRHWS calculations, the neutron rms radius 
ry and Fermi energy An of 72Ca converge around Rmax 
= 35 fm while a independence of the binding energy per 
nucleon £/A and proton rms radius rp on the box size 
can be seen. These sets of parameters in SRHWS, Ecut 
= 75 MeV and Rmax = 20, 25, 30 and 35 fm, correspond 
to cutoff’s on principal quantum number Nmax = 25, 31, 
37 and 43 which are used in SRHHO calculations in order 
to make fair comparisons between SRHWS and SRHHO. 
Similar as those from SRHR and SRHWS, £/A and rp 
depend little on Nmax in SRHHO. However, the neutron 
rms radius rn increases with Nmax steadily, which shows 
a much slower convergence. As it is based on a complete 
basis, SRHHO can reach convergence of ry if Nmax is 
large enough. From Table V, one finds that for the same 
Nmax (or equivalent Rmax), a difference of Ar, ~ 0.2 
fm between SRHHO and SRHWS (SRHR) can be seen. 
From the slow convergence of rn with Nmax in SRHHO 
(ANmax = 6 gives Ar, ~ 0.02 fm), we can estimate the 
lower limit of the Nmax aS Nmax œ 90 in order to give ry 
= 4.8 fm in SRHR or SRHWS. 
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FIG. 9: Comparison of density distributions for Ca from 


SRHR, SRHDWS and SRHHO. The parameter set NLSH 
is used for the Lagrangian. Ar = 0.1 fm for SRHR and 
SRHDWS. Ecut = 75 MeV and Rmax = 20 (thick dot-dashed 
curve), 25 (thick long-dashed curve), 30 (thick dashed curve) 
and 35 fm (thick dotted curve) for SRHDWS. In SRHHO 
calculations, Nmax = 25 (thin dot-dashed curve), 31 (thin 
long-dashed curve), 37 (thin dashed curve) and 43 (thin dot- 
ted curve) are used. The density distribution from SRHR are 
very close to that from SRHSWS with the same box size. For 
brevity, only pn(r) from SRHR with Rmax = 35 fm (thick 
solid line) is displayed here. 


We compare the neutron density distribution of 7?Ca 
from different SRH approaches in Figs. 8 and 9. With 
the same box size, the density distribution from SRHR 
are almost identical with those from SRHWS, which in- 
dicates the equivalence between SRHWS and SRHR. For 
brevity, only pa(r) from SRHR with Rmax = 35 fm is 
displayed in Figs. 8 or 9 which covers the curve corre- 
sponding to p,(r) from SRHSWS with Rmax = 35 fm 
in Fig. 8. On the other hand, pn(r) from SRHHO even 
with Nmax = 43 fails to reproduce the result of SRHR 
due to the localization property of HO potential. In addi- 
tion, with the same Nmax, the spatial extension of pa (r) 
from SRHWS are always larger than that from SRHHO. 
The variational tendency of the curve py(r) ~r also ex- 
plains different convergence behaviors of ry, in SRHWS 
and SRHHO as given in Table V. With increasing Rmax, 
the py(r) in SRHWS has the correct asymptotic behavior 
while that in SRHHO decay too quickly. 


This result is very encouraging and tell us that even 
the long tail (or halo) behavior in neutron density dis- 
tribution for nuclei near the drip line can be described 
by SRHSWS as well as SRHR, if pairing correlations is 
incorporated properly. 
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TABLE VI: Effects of negative levels on bulk properties in 
SRHDWS for '®O. The parameter set for the Lagrangian is 
NLSH, Rmax = 20 fm, Ar = 0.1 fm, and Nj, = 25. For the 
initial Woods-Saxon like potentials, parameters in Ref. [23] 
are used except for Vo which is specified in the table. The left 
value in each entry gives the result without negative levels 
included and the right one that with Nmax = 25. Energy is 
in MeV and radius in fm. 


Vo = 54 MeV Vo = 72 MeV Vo = 90 MeV 
EJA 8.547 | 8.013 8.117 | 8.015 8.427 | 8.012 
Trms 2.385 | 2.568 2.531 | 2.567 2.610 | 2.567 


C. Contribution from negative levels in the 
SRHDWS theory 


In the expansion of the nucleon wave function, 
Eq. (21), one has to take into account not only the levels 
in the Fermi sea but also those in the Dirac sea because 
they form a complete basis together. We study the ef- 
fects of negative levels and give results of an example, 
160, in Table VI. Firstly, without negative levels in- 
cluded, the nucleus is over bound and the nuclear size 
is smaller as seen from Table VI. Secondly, in contrary 
to the case with negative levels included, the calculated 
nuclear properties depend on the initial potentials very 
much if no negative levels included. 


It should be noted that the contribution from negative 
levels depends on the initial Woods-Saxon potentials for 
generating the DWS basis. So do the cutoff Nmax or Ecut 
for convergence. If the initial Woods-Saxon potential is 
exactly identical to the converged potentials, the matrix 
in Eq. (22) is diagonal, negative states do not contribute 
because of the no sea approximation. Positive states can 
also be chosen as less as possible, e.g., 151/2, 1p3/2 and 
1pı/2 are enough for 16O. From the third column corre- 
sponding to Vo = 72 MeV in Table VI, one finds that the 
initial nuclear potential for the Dirac equation proposed 
in Ref. [23] is a good choice for SRHDWS as the negative 
states only contribute ~ 1.25 % to both E/A and rims. 
If we change the initial potentials, e.g., by changing Vo 
by 25%, much larger contributions from negative states 
are found in Table VI. 


In order to know the contribution of negative levels in 
the Dirac sea to the wave function, the value of `, |c} | 
in the expansion, Eq. (21) has been given in Table VII for 
occupied states of 160. We note that the nucleon wave 
function is normalized to one. It can be seen that a small 
component from negative states in the wave functions, 
about 1074~75, contributes to the physical observables 
such as £/A and ryms by the magnitude of 1%~10% as 
seen from Table VI. Again we notice that the initial 
Woods-Saxon potentials differ more from the converged 
ones, the larger is the contribution from negative levels. 


TABLE VII: Contribution of negative levels in the Dirac sea 
to the single nucleon wave functions of ‘©O in SRHDWS. For 
each single particle levels, X |c} |? in the expansion Eq. (21) 
is presented. The parameter set for the Lagrangian is NLSH, 
Rmax = 20 fm and Ar = 0.1 fm and Noax = Nmax = 25. 
For the initial Woods-Saxon like potentials, the parameters 
in Ref. [23] are used except for Vo which is specified in the 
table. The binding energy per nucleon and nuclear radius for 
each calculation is presented in Table VI. 


Level Vo = 54 MeV Vo=72MeV W=90 MeV 
Nisız 915x10” 348x10” 344x107 
N1p3/2 5.32x 107° 6.05x 107° 3.69x 10-4 
N1pi/2 4.01x10~4 1.06x10~° 8.80x10~4 
Pls1/2 1.92x 10-4 4.92x107° 4.89x107° 
Plps/2 1.17x10~4 3.27x107° 8.07x10~° 
P1p1/2 8.09x 10-4 2.08 1074 2.38x107° 


V. SUMMARY 


We have solved the spherical relativistic Hartree the- 
ory in the Woods-Saxon basis (SRHWS). The Woods- 
Saxon basis is obtained by solving either the Schrodinger 
equation (SRHSWS) or the Dirac equation (SRHDWS). 
Formalism and numerical details for both cases are pre- 
sented. The WS basis in the SRHDWS theory could be 
very smaller than that in the SRHSWS theory. This will 
largely facilitate solving the deformed problem. 

The results from SRHWS are compared with those 
from solving the spherical relativistic Hartree theory in 
the harmonic oscillator basis, SRHHO, and those in the 
coordinate space, SRHR. For stable nuclei, all approaches 
give identical results for properties such as total binding 
energies and the neutron, proton and charge rms radii as 
well as neutron density distributions. 

For neutron drip line nuclei, e.g. Ca, which has a 
very small neutron Fermi energy A, ~ 0.2 MeV, both 
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SRHR and SRHWS easily approach convergence by in- 
creasing box size while SRHHO does not. Furthermore, 
SRHWS can satisfactorily reproduce the neutron density 
distribution from SRHR, but SRHHO fails with similar 
cutoff’s. 

In SRHDWS calculations, the negative levels in the 
Dirac Sea must be included in the basis in terms of 
which nucleon wave functions are expanded. We stud- 
ied in detail the effects and contributions of negative 
states. Without negative levels included, the calculated 
nuclear properties depend on the initial potentials very 
much. A small component from negative states in the 
wave functions, about 10~*~~°, contributes to the phys- 
ical observables such as Æ/A and ryms by the magnitude 
of 1%~10%. When the initial potentials differ more from 
the converged ones, the contribution from negative levels 
becomes more important. 

We conclude that the Woods-Saxon basis provides a 
reconciler between the harmonic oscillator basis and the 
coordinate space which may be used to describe exotic 
nuclei both properly and efficiently. 

The extension of spherical relativistic Hartree theory 
in the Woods-Saxon basis to deformed cases with pairing 
included is in progress. 
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